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We show that a dissipative current component is present in the dynamics generated by a Liouville- 
master equation, in addition to the usual component associated with Hamiltonian evolution. The 
dissipative component originates from coarse graining in time, implicit in a master equation, and 
needs to be included to preserve current continuity. We derive an explicit expression for the dissipa- 
tive current in the context of the Markov approximation. Finally, we illustrate our approach with a 
simple numerical example, in which a quantum particle is coupled to a harmonic phonon bath and 
dissipation is described by the Pauli master equation. 
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Open quantum systems are encountered in many do- 
mains of physics, ranging from condensed matter physics 
0, 0] to quantum optics 0, chemical physics Q and 
NMR theory In all these fields one is interested in 
the reduced dynamics of a quantum system in contact 
with its environment. The reduced dynamics is usually 
approximated by a master equation. A simple example 
is the Pauli master equation (PME) [fj in which the di- 
agonal elements Su(t) = (J, \S(t) \ l) of a coarse grained 
reduced density operator S(t) in the basis of the system 
eigenstates evolve in time according to: 

7 Q 

-jf = J2 (WlmSmrn ~ W ml S U ) . (1) 

The transition probabilities Wij are given by Fermi's 
golden rule expression in terms of the coupling to the 
bath. The Wij satisfy detailed balance, which guaran- 
tees that the thermal equilibrium distribution is the sta- 
tionary solution of Eq. QJ. A difficulty with Eq. is 
that it describes a charge (or probability) density evolu- 
tion that does not seem to satisfy the current continuity 
equation, given that the current associated with eigen- 
states is divergenceless. This apparent violation of cur- 
rent continuity makes the use of Eq. JTJ problematic to 
model electron transport kinetics [7J, |8( , although it has 
been argued that under steady state conditions Eq. (JJJ 
can still be used to simulate electron transport in very 
small devices Q. 

Eq. is a special case of a Liouville master equation 
of the form: 

^. = -i[H,S(t)]+C[S(t)l (2) 

where H is the Hamiltonian of the system and C[S(t)] 
describes the dissipative part of the dynamics after elim- 
inating bath memory effects by coarse-graining in time 



(Markov approximation) 0,^3- Eq. © is derived by as- 
suming a harmonic bath and by treating the system-bath 
interaction to 2nd order of perturbation theory. The lim- 
its of validity of such approximation and of the Markov 
approximation requiring separation of the time scales of 
system and bath, have been amply discussed in the liter- 
ature IelUI- 

Eq. (J2J) applies to a generic interacting N— particle 
quantum system described by a N— particle coarse- 
grained reduced density operator S. 

In what follows we adopt atomic units (e = m = 
h = 1). We indicate coarse-grained quantities by a 
bar. The coarse-grained charge density n(x; t) is given 
by n(x; t) = Tr {S(t)p(x)} in terms of the charge density 
operator p(x) = J2i=i n^( x — ^0- The coarse-grained 
Hamiltonian current density is given by jH(x;t) = 

Tr (B(t) J (x)\, in terms of the current density opera- 
tor j(x) = \ Y,i=i,n [PiH x - &i) + 5(x - xi)pi]. If the 
potential in H is local in space, then the continuity equa- 
tion n(x\ t) = — V • Jh{x\ t) is evidently satisfied if 7]: 

Tr{C[S(t)}p(x)}=0. (3) 

Generally master equations in which the dissipative dy- 
namics leads asymptotically to the thermal equilibrium 
distribution do not satisfy Eq. © lj| • This is at vari- 
ance with the exact reduced dynamics which satisfies the 
continuity equation if the coupling to the bath is local. 

In this letter we show that assertions on the violation of 
continuity by master equations that do not satisfy Eq. (J3J) 
are indeed incorrect, and that the coarse-grained (Marko- 
vian) dynamics exactly satisfies current continuity, as the 
true reduced dynamics does if the coupling to the bath 
is local. This is because the following equation holds: 

T*{C[S(t)]p(x)}=fic(x) = -V-Jc(x;t) (4) 
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In Eq. jo is a dissipative current originating from the 
action of the bath on the system on the coarse grained 
time scale of Markovian dynamics. This current contri- 
bution, which has been ignored so far, is of the same order 
of the Hamiltonian current in strongly inhomogeneous 
systems. We also provide an expression for jc, whose 
calculation only requires knowledge of the coarse-grained 
density operator S(t), of the eigenvalues and eigenvectors 
of the system Hamiltonian H , of the Hamiltonian R of 
the harmonic reservoir, and of the coupling K between 
the system and the bath. 

The derivation of a master equation of the form @ 
can be found in textbooks [lOj. Here we limit ourselves 
to pointing out the key elements of this derivation that 
are important for our discussion of current continuity. 

Let S(t) be the reduced density operator in ab- 
sence of coarse graining. As usual, to separate the 
effect of the bath from system dynamics, it is conve- 
nient to work in the interaction picture (IP): S(t) — 
e -tH(t-t ) s (£j e iH(t-t ) ^ wnere f Q j s a time preceding t and 

s(t) is the reduced density operator in the IP. Taking the 
time derivative of S(t), the reduced density operator in 
the Schrodinger picture (SP), one gets: 



dS{t) 
dt 



-i[H,S(t)} 



Eq. JSJ) governs the exact dynamics of the system. 
The first term on the right hand side (r.h.s.) corre- 
sponds to the usual Hamiltonian dynamics, which sat- 
isfies current continuity if the potentials in H are lo- 
cal. For the second term on the r.h.s. one can show that 
Tr je-^-'o) ilMl e iH{t-t°) p(x)\ = 0, if the interaction 

potential with the bath is local. Therefore, this term does 
not modify the instantaneous distribution of charges, and 
the exact dynamics satisfies current continuity. 

In order to convert Eq. (0) into a master equation of 

in 

(to) 



the form J5J, one replaces the time-derivative 



s (t)- 



by a finite difference approximation - 
calculated to leading order of perturbation theory using 
Fermi's golden rule, and then one takes the limit At — > 
(corresponding to coarse- graining in time) in the equa- 
tion for system dynamics, i.e. one writes for the coarse- 
grained density operator S(t): 

= -* \H,S(tj\ + lim e -mt-t )^ e iH(t-t ) 
dt to->t At 



[H,S(t)} + lim 



As 



t ->t At 



(6) 



We recall that Fermi's golden rule expression for A| 
is independent of At. The limit for At — > is possible 
because At, although large relative to inverse system fre- 
quencies for Fermi's golden rule to be applicable, is small 
compared to the relaxation time. Thereafter all reference 
to to disappears from Eq. @ and the future evolution of 
the system is entirely determined by its present and not 



by its past (Markov approximation) . Eq. © is the same 
as Eq. © 

To prove that Eq. (J2J satisfies charge continuity even 
when Eq. y| is not satisfied, we compute the time deriva- 
tive of the coarse-grained charge density: 



dt I dt 



(J) 



= -zTr { [H, S(t)] p{x)} + Tr /lim ^P(*)} ■ 

The first term on the r.h.s. of Eq. @ gives: 
-iTi{[H,S(t)] p(x)} = -V ■ J H (x;t). To compute 
the remaining term we consider the exact density op- 
erator S(t') = e~ iH{ ? '- to h(t')e lH(r '~* o) , whose dynam- 
ics includes the effect of the bath, and the density op- 
erator Suit'), whose dynamics is governed by H alone 
and does not include the effect of the bath: Sji(t') = 
e -iH(t -t ) s ^ e iH(t -to)^ Here we have assumed that at 
time t S(t) — Su(t)- Then at time to, S(t ) = s(t ) and 
S H (to) = s(t). Thus: 



As = s(t) - s(t ) = S(t) - SQto) - S H (t) + S H (t ) 



dt 



la 



' (s(t')-s H (t': 



(8) 



Since both the exact and the H— dynamics satisfy conti- 
nuity, we can write 



Tt{Asp(x)} 



dt' (n(x;t') - n H (x;t')) 



in 



-V- / dt' (j(x;t')-j H (x;t')) 
'to 
rt 



-V- f dt'j c (x;t'), 

J to 



(9) 



where j H = Tr {Sjjj} and j = Tr [sj}- In Eq. © the 
difference between the exact current, which includes the 
effect of the bath, and the if— current defines the current 
jc = j — jH- Using © we obtain the second term on the 
r.h.s. of Eq. Q: 

= - v -(fes*jC*' ic(x;t) )- (io) 

We identify the term in parenthesis in Eq. I|1U|) with the 
coarse-grained current jc. Thus we obtain 



dn(x) 
dt 



-V-j H (x;t)-V-jc(x;t) = -V-j{x;t). (11) 



Eq. JTTJ is the main result of this Letter: it shows 
that the master equation does not violate current conti- 
nuity and that the coarse-grained physical current J, i.e. 
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the current measured in experiments, is the sum of two 
contributions, a current Jh , associated with the charge 
flow from Hamiltonian propagation in the Liouville mas- 
ter equation J2J, and a current jc, associated with the 
charge flow due to inelastic collisions with the bath in 
the same equation. In the following we show that the 
limit in Eq. (|10|l can be calculated to leading order of 
perturbation theory, i. e. at the same level of approxima- 
tion used in the derivation of the master equation itself. 
The resulting expression for Jc is independent of to- Thus 
j(x; t) can be calculated from the knowledge of S(t) at 
time t only. 

For the derivation of an explicit expression for jc , we 
consider a coupling K of the form K = J2 a V a F a , where 
V a = J2i=i n Wi^i) i s a l° ca l potential acting on the 
system, and F a acts on bath variables only. The sum 
is over all the bath eigenmodes (phonons), but to sim- 
plify the notation we drop the index a in what follows. 
A standard procedure [lfj gives As to second order of 
perturbation theory: 



As = s(t) - s(t ) 



B 



mn,opW mn 



ninop 

rt-to 



d£e 



(12) 



where the sums are over eigenstates \i) of H having en- 
ergy ej. The operators A and B are: 



(n \o) V op (jp\ s(t ) 
- \o) V op (p\ s(t ) \m) V mn (n\ 

Bmn.op = \o) V op (p\ s(t ) \m) V mn H 

-s(t ) \m) V mn (n |o) V op (p\ , 
where V mn — (m \V\ n), and 



(13) 



e i(e m -e n )f (F(t')F(Q)) R dt', 
e i(e -e„)t' (i?(0) J F(t')) i? dt' 



(14) 



are bath spectral densities, which, for a harmonic bath, 
can be calculated analytically and are zero when m = n 
(o = p). In Eq. lO F(t') is in the IP and (• • -) R denotes 
an average over the bath degrees of freedom. If t — to is 
large compared to the inverse level spacing 1/ (ej — ej), the 
time integral in Eq. 1|12J) leads to selection rules, where 
the only terms that survive are those with m = p and n = 
o . The resulting As varies linearly with At . Inserting 
■^f in Eq. ©, one obtains the usual expression for the 
dissipative master equation |10| . 

In order to compute the dissipative current, one has to 
evaluate 



j c (x,t)= lim -J- / dt'Tv{(S(t')-S H (t'))J(x)} 



Using (S(f) -S H (t')) 

e -iH(t'-to) ( s (^) _ s ( tQ ) _ ( s (i) _ s ( <0 ))) e iK(t'-t ) 5 

one can insert the expression l|12l) for s(t') — s(t ) and 
s(t) — s(t ), respectively, into Eq. I|15|). 

The additional time-integral and the exponentials lead 
to different selection rules than in case of Eq. l|12ll . Iden- 
tifying the leading terms, and finally letting to — » t , we 
obtain an expression for the dissipative current: 



1 



{T kl -Tl k i)J{x)i k , (16) 



Jc{x;t) = i 



where J{x)i k = (l\J(x)\k) and U k u%i are given by: 



fckl = - ^2 \Vnk\ 2 S k l(t)w nk 

-J2s kl (t)\V ln \ 2 w+, 
%i = - 1 ^V kn V n iSu(t)w~ l 

+ 2 V kn S nn V n i(w+ + w^J 

n^k 

— 53 s kk (t)v kn v r , 



(17) 



i w L 



(18) 



(15) 



Eq. (|16|l is the second important result of this pa- 
per and completes the proof that jc(x;t) can be calcu- 
lated from knowledge of S(t) and microscopic Hamilto- 
nian only. Eq. l|16f) gives (to 2nd order of perturbation 
theory) the coarse-grained current associated to phonon 
induced inelastic transitions between all pairs of system 
states. The transition matrix elements depend on (am- 
plitude and phase of) the coupling potential and on the 
fluctuations of the bath (via the bath spectral densities). 
The derivation of Eq. I |16(l does not require additional 
assumptions besides those made in the derivation of the 
Liouville master equation. In this context Eq. (|16J) is 
exact. Our approach requires an explicit phonon (bath) 
model, however simplified. 

Sometimes master equations are used without refer- 
ence to an explicit bath model. This is the case, for in- 
stance, with the common relaxation time approximation. 
This approximation does not violate charge continuity 
when it is applied to spatially homogeneous systems. In 
inhomogeneous systems, however, the relaxation time ap- 
proximation leads to violation of charge continuity when- 
ever the bath originates a spatially non-uniform instan- 
taneous change of the coarse-grained charge density, i.e. 
whenever hc(x) ^ [Tl) . 

Finally, we use a simple numerical example to illustrate 
our theory. In this example a single quantum particle 
(electron) is in a one dimensional potential U(x). The 
electron is initially placed in the first excited eigenstate 



4 




Position 



FIG. 1: Upper Panel: Calculated charge densities of the 
ground state (solid line) and first excited state (dashed line). 
Middle Panel: Dissipative current Jc(x) at a given instant 
of the propagation, calculated from Eq. I|16|l . Lower panel: 
Divergence of jc(x) (solid line), time derivative of the charge 
density n(x; t) at the same time (dashed line), and the sum of 
the two quantities (dotted line). It is evident that continuity 
equation is satisfied during relaxation. 



of the Hamiltonian H = — ^ + U (x) . In absence of 
dissipation the electron would stay in the excited state, 
but due to coupling with a phonon bath it relaxes toward 
thermal equilibrium. In this case, neglecting the imagi- 
nary part of the bath spectral densities, Eq. J2J) reduces 
to Eq. PJl. i.e. the PME, and the density matrix stays 
diagonal while the system approaches equilibrium. Even 
though the density matrix does not carry any current, 
Eq. 1|16|) gives a non-zero current due to the Tj~i terms 
during relaxation dynamics. 
We take 



w ' = u> 



n(e r , 
n(e r . 



- e r . 



> e Tl 
< e,, 



(19) 



where n(cu) = l/(efcr — 1) is the mean occupation num- 
ber of phonons with energy uj. This simple choice guar- 
antees detailed balance. We assume that V = dv } x ' is 

ax 

the potential that characterizes the coupling with the 
bath. The transition probabilities in Eq. Q are given 
by Wi m = 2\Vi m \ w^. We represent the space coor- 
dinate x on a fine numerical grid with 256 points. In 
Fig. ^ we report the calculated charge density n Q (x) of 
the ground state, the calculated charge density ni(x) of 
the first excited state, and the calculated time derivative 
of the charge density n(x] t) at a time t during relaxation 
dynamics (in this example we take kT = 0.4 (ei — eo)). 
It is clear that there can be no current associated to the 
eigenstates to compensate the change of charge density 
n(x] t). 

In the same figure, we report the dissipative current 



density jc(x;t) at the same time t. It is obtained from 
Eq. (|TB|> . where TZki = because S(t) is diagonal. As we 
can see, the continuity equation n(x;t) + -^Jc(x;t) = 
is exactly satisfied within the numerical accuracy of the 
calculation. 

Numerical integration of Eq. © for a single quantum 
particle is feasible, but is not possible in general for in- 
teracting many-body quantum systems. In practical ap- 
plications to a many-body quantum system, the latter 
is usually approximated by a system of non-interacting 
particles in an effective single-particle Hamiltonian like 
e.g. the time dependent Hartree Hamiltonian Then 
Eq. J5J, the N— particle Liouville master equation, can be 
reduced to an equation for an effective single-particle den- 
sity operator in which the master term depends quadrat- 
ically rather than linearly on the (single-particle) den- 
sity operator [J^. The basic assumptions on perturbative 
coupling to the bath and coarse-graining in time that led 
to Eq. 

remain valid, however, and so is our analysis 
that led to Eq. <|Td|| . Also in this case, a proper expression 
for jc(x;t) can be found, which restores exactly current 
continuity removing a major limitation to the ap- 
plication of quantum kinetic approaches to the study of 
electron transport problems. 
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